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Abstract 

To investigate non-linear dynamical systems, like for instance artificial 
satellites, Solar System, exoplanets or galactic models, it is necessary to 
have at hand several tools, such as a reliable dynamical indicator. 

The aim of the present work is to test a relatively new fast indicator, 
the Mean Exponential Growth factor of Nearby Orbits (MEGNO), since it 
is becoming a widespread technique for the study of Hamiltonian systems, 
particularly in the field of dynamical astronomy and astrodynamics, as well 
as molecular dynamics. 

In order to perform this test we make a detailed numerical and statistical 
study of a sample of orbits in a triaxial galactic system, whose dynamics 
was investigated by means of the computation of the Finite Time Lyapunov 
Characteristic Numbers (FT-LCNs) by other authors. 
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1. Introduction 

In the present work we accomplish an exhaustive study of the MEGNO 
when applied to a given sample of orbits in a triaxial galactic potential studied 
by [H- In that work, the authors use a well-known tool, the Lyapunov 
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Characteristic Numbers (see e.g. |2j), to identify the character of the selected 
orbits in order to classify them as regular or chaotic. 

The MEGNO is introduced by ^ and, in jij, this technique is formalized 
and its application extended to discrete Hamiltonian systems like maps; also 
a generalization of the MEGNO is introduced therein. This tool has become 
of widespread use for studying several astronomical problems as well as many 
other Hamiltonian systems (see, e.g. 
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the MEGNO succeed in furnishing a clear insight of the 
global structure of the phase space of simple multidimensional Hamiltonian 
systems, providing a clear picture of the resonant network as well as the 
regular and chaotic domains. 

Herein instead, a far more complex non-linear system is addressed that 
reproduces many characteristics of real elliptical galaxies, namely, the one 
introduced by |l|. This model will be used as the scenario for a detailed 
comparison between the MEGNO and the Lyapunov Characteristic Numbers 
and even the Fast Lyapunov Indicator (FLI) introduced by (20[ . 

There are many efficient dynamical indicators, some of them based on de- 
viation vector (s), for instance, the spectra of stretching numbers, helicity and 
twist angles, the computation of the alignment indices introduced by Skokos, 
the Relative finite time Lyapunov Indicator (RLI) and also the Avera ge P ower 
Law Exponent (APLE), a technique recently developed 2l|, |22|, |23j, |24|, and 
others based on spectral analysis, such as j25|, the Frecuency Map Anal- 
the one due to Sidlichovsky and Nesvorny 
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[28| . and perhaps 

the latest one, the FMI (Frequency modulation indicator) [29]. However the 
present paper is devoted to accomplish a thorough test of the MEGNO, point- 
ing out not only its advantages but its drawbacks as well. Therefore herein we 
just focus our attention on an exhaustive comparison of the MEGNO against 
the Lyapunov Characteristic Numbers (and eventualy the FLI), since with- 
out any doubt, the latter is the most widespread tool in, at least, the last 
forty years, and it is still being used by many authors. On the other hand, 
as far as we know, a full test of the MEGNO when applied to a non-linear 
somewhat realistic Hamiltonian system has not been performed yet. 



2. The Mean Exponential Growth factor of Nearby Orbits (MEGNO) 

In this section we summarize the main features of the MEGNO (described 
in detail in jij). This is an alternative tool to explore the phase space which 
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belongs to the class of the so-called fast indicators. 

Let H (p, q) with p, q G M, N be an iV-dimensional Hamiltonian, that we 
suppose autonomous just for the sake of simplicity. Introducing the following 
notation: 

x = (p, q) G R 2N , v = (-dH/dq, OH /dp) G R 2N , 
the equations of motion can be written in a simple way like 

x = v(x). (1) 

Let 7(x ; t) be an arc of an orbit of the flow ([T]) over a compact energy 
surface: C M? N , = {x : H(p, q) = h} with h = constant, then 

7 (x ; t) = {x(t' ; xo) : x G M h , < t' < t}. 

We can gain fundamental information about the Hamiltonian flow in the 
neighborhood of any orbit 7 through the largest Lyapunov Characteristic 
Number (LCN) defined as: 

(7(7) = lima 1 (7(x ;t)), (7 1 ( 7 (x ; t)) = - In [||5 7 (x ; t) \\] , (2) 

with 57(x ; t) an "infinitesimal displacement" from 7 at time t, where || • || is 
some norm. The fact that the LCN measures the mean exponential rate of 
divergence of nearby orbits it is clearly understood when Eq. (j2J) is written 
in an integral fashion: 

1 fhfarf) 



a(j) = lim - / 7 f. dt' =(6/6), (3) 

where 6 = \\6\\, 6 = 6.5/ dt = 5 ■ 5/\\5\\, and the bar denotes time average. 
Also, the tangent vector 5 satisfies the variational equation 

5 = A( 7 (x ; £))■<?, 

where A is the Jacobian matrix associated with the vector field v. 

Now we are in a position to introduce the MEGNO, y(7(x ; t)), through 
the expression: 

2 f'hi^t').^ 



2 f H^tT) 
y(7(x ;t)) = - / , A 6t 
t Jo <57(x ;t') 
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which is related with the integral in Eq. ([3]); i.e., in case of an exponential 
increase of 5, <57(x ; t) = ^7(x ; to) ■ exp(At), the quantity F(7(x ; t)) can be 
considered as a weighted variant of the integral in Eq. fl3]) . Instead of using 
the instantaneous rate of increase, A, we average the logarithm of the growth 
factor, In [57(x ; t)/Sj(xo] t )] = Xt. 

Let us describe the MEGNO's asymptotic behavior to exhibit its ability 
to give a clean idea of the character of orbits. Firstly, consider the case of 
orbits on irrational tori for a non-isochronous system. As it is shown in 
for quasi-periodic orbits, j g , the temporal evolution of Y( n f q (xo;t)) is given 
by ' " 

Y ( 7g (x ; t)) « 2 - ln(1 ^ gt)2 + O ( 79 (x ; t)) , 

A q I 

where X q is the linear rate of divergence around j q and O is a null average 
oscillating term. Accordingly to this formula, the lirn^oo Y (7 ? (x ; t)) does 
not exist, but on introducing a time average 

t 



Y{ lq {^ 0] t)) = - t J y^xo;*'))^' 



it can be found that 

Y{ lq )= limF(7,(xo;f)) = 2. 



t— >oo 



Then, for quasi-periodic motion, Y(y) is a fixed constant, independent of 7. 

When taking irregular orbits 7«, i.e. orbits on some stochastic layer, for 
which <57i(x ; t) ~ 57(x ; to) • exp (ait), ai being the LCN of 7*, the temporal 
evolution of the MEGNO will be given by: 

Y ( 7 i(x ; t)) ^a,t + ( 7l (x ; t)) , 

with O some bounded amplitude and null average oscillating term (see 0]). 
On averaging over a sufficiently large interval we have: 

F(7i(x ;t)) w y t, t->oo. 

Therefore, in the case of chaotic orbits, not only Y (7i(x ; t)) but also Y (71 (x ; t)) 
grow linearly with time, with a slope equal to the LCN of the orbit or one 
half of it, respectively. Wherever the phase space has a hyperbolic structure, 
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Y will indefinitely grow with time. Otherwise, it will approach a constant 
value, even in the degenerated case in which S grows with some power of t, 
e.g. n, for which Y — > 2n when t — > oo. 

We notice that the temporal evolution of the MEGNO can be briefly de- 
scribed in a suitable and unique expression for all kind of motion. Indeed, the 
asymptotic behavior of F(7(x ; £)) can be summarized in the following way: 
F(7(x ; £)) ~ a 7 i + g? 7 , where a 7 = cr 7 /2 and d 7 ~ for irregular, stochastic 
motion, while a 7 = and d 7 ~ 2 for quasi-periodic motion. Deviations from 
the value cf 7 « 2 indicate that 7 is close to some particular objects in phase 
space, being d 7 < 2 or d 7 > 2 for stable periodic orbits (or resonant elliptic 
tori), or unstable periodic orbits (or hyperbolic tori) respectively (see [i| for 
details). Finally, the quantity &i =Y/t verifies that 

2 

(7i(7g(x ;t)) « -, 0-1(7^X0;*)) « o-j, as t -)■ 00, 

supporting the fact that in regular domains, &i converges to faster than a x 
(which goes to zero like lnt/t), while for stochastic domains, both quantities 
tend to the positive LCN at a rather similar rate. 

Let us introduce here a brief comment regarding the computation of the 
LCN. As it is already well-known, though the definition of the Lyapunov 
Characteristic Numbers encompasses an integration over an infinite interval 
of time, their numerical computation involve a rather large but finite time 
interval and the expected null value corresponding to regular motion is un- 
likely to be reached. In such instead, the 'Finite Time Lyapunov 
Characteristic Numbers' (FT-LCNs hereafter, following the nomenclature 
given in |30|) attains a value of order In T/T, being T the total integration 
time. Thus, a critical value has to be adopted as 'zero', so that FT-LCNs' 
values greater or lower than such critical value are regarded as different from 
or equal to zero respectively. 

3. The potential 

For the comparison of the MEGNO vs. the FT-LCNs we deal with the 
potential introduced by [l[ which, obtained after the virialization of an iV- 
body self-consistent model composed of one hundred thousand particles, re- 
produces many features of real elliptical galaxies, such as mass distribution, 
flattening, triaxiality and rotation (see also, |3l|). Nonetheless, it is clear 
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that a real elliptical galaxy is a much more complex astrophysical system 
than a purely dynamical one. 

This potential seems to provide an adequate scenario for the compar- 
ison between the two above mentioned techniques. To this aim, we ad- 
dress the study of the set of randomly selected orbits O = {xj(t),i = 
1, ■ • • , 3472, Xj(0) = x j}, classified by means of the FT-LCNs in [lj. Their 
initial conditions x j and their concomitant FT-LCNs values were provided 
by the authors. We will identify each orbit of the set O with a label that 
runs from 1 to 3472. 

The equation that reproduces the potential is: 



V(x, y, z) = -f (x, y, z) - f x (x, y, z) ■ (x 2 - y 2 ) - f z {x, y, z) ■ {z 2 - y 2 ), (4) 
where 

fn{ X lVl Z ) = ~ " Tac5~; (5) 

\Pn n + d n"\ an 



where "p 2 n is the square of the softened radius given by p 2 l = x 2 + y 2 + z 2 + e 2 
when n = 0, or p 2 n = x 2 + y 2 + z 2 + 2 • e 2 for n — x, z, and a n , S n , a n , ac n 
are constants. The adopted value for the softening parameter is e ~ 0.01 
for any n. The functions f n (x,y,z) were computed through a quadrupolar 
iV-body code for a hundred thousand bodies, which allowed the authors to 
write them in a general fashion given by Eq. (JSJ) . The adopted values for the 
constants a n , 5 n , a n and ac n are given in Table [T] (further references in [32]). 
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a 


5 


ac 


n 


= 


0.92012657 


1.15 


0.1340 


1.03766579 


n 


= X 


0.08526504 


0.97 


0.1283 


4.61571581 


n 


= z 


-0.05871011 


1.05 


0.1239 


4.42030943 



Table 1: Adopted values for the coefficients of the functions f n given by Eq. ([5]). 

The stationary character of the parameters given in Table [1] were tested 
by performing several fits at different times after virialization, resulting with 
a precision of 0.1%. 
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After the system had relaxed, there remained 86.818 bodies resembling 
an elliptical galaxy (the system obeying a de Vaucouleurs law, as Fig. 2 in 
[l[ shows) with a strong triaxiality and a flattening that increases from the 
border of the system to its center (see Table I in the same work). 



Fig. [U taken from 32J, displays the behavior of the f n regards to r, being 
fz < 0, /o, fx > 0, and f x > \f z \ for the whole r range, while for r > 0.36 
it is fo > fx- Notice that the functions f n are plotted with the concomitant 
sign with which they appear in Eq. (JJJ). The obtained triaxial potential has 
semi-axis X, Y, Z satisfying the condition X > Y > Z and its minimum, 
which is close to —7, matches the origin. The potential is less flattened than 
the mass distribution, as expected (see Table I in [l[). 



-1G 



Figure 1: Radial dependency of the functions /„ with their concomitant signs (figure taken 
from [33 |. 



4. Comparison of the MEGNO vs. the FT-LCNs 

The present section is devoted to performing a numerical and statistical 
comparative study of the results obtained by recourse to the MEGNO when 
applied to the set O of orbits and those provided by [lj through the FT-LCNs. 

It is of interest to deem the energy spectrum of the orbits in the set O 
displayed in Fig. |2} where we observe that most of the 3472 orbits considered 
have large energies; indeed, in the main they have energies in the range 
-3 < E < 0. 

Let us recall the criterion used in [l] to classify the orbits in O according to 
their FT-LCNs: those orbits with their largest FT-LCN below some critical 
value V c were labeled as regular, otherwise they were classified as chaotic. 
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Figure 2: Energy spectrum of the 3472 orbits. 



Since the numerical integrations carried out by |l| for the computation of 
the FT-LCNs encompasses an interval of 10000 u.t. (units of time), the ex- 
pected value for V c would be V* = In T/T « 0.00092 (u.t.) -1 . Notwithstand- 
ing, [l[ took an empirical value slightly higher, = 0.00155 (u.t.) -1 , and this 
is the one we consider in order to observe their classification into regular and 
irregular orbits. Along this investigation we adopt for the MEGNO a thresh- 
old value of 2.01 for regular orbits. According to pj, the set O = Oc U Or, 
where Oc and Or include 1828 chaotic orbits and 1644 regular ones, respec- 
tively. 

The computation of the MEGNO, as well as that of the largest FT-LCN, 
requires the integration of the equations of motion along with their first 
variationals, the initial conditions for the latter being taken at random in 
phase space and with unit norm. The integrations were accomplished using 
a R ung e-Kutta 7/8 th order integrator (the so-called DOPRI8 routine -see 
33ll34j-). over short: 5000 u.t, intermediate: 10000 u.t. and large integration 
times: 100000 u.t. The precision in the conservation of the energy was of the 
order of ~ 10~ 12 . 

The FT-LCN values corresponding to a total integration time of 10000 
u.t. for the sample of orbits classified in [lj were kindly provided by Muzzio, 
to whom we are grateful. 

Herein we present the results corresponding to 5000 u.t. and to 10000 
u.t. and even larger motion times, in order to disclose how efficient could 
this tool be, to provide dynamical information at short times. 
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Figure 3: Correlations between 2Y/T and the largest FT-LCN for chaotic orbits (left 
panel) and between Y/T and the largest FT-LCN for regular orbits (right panel) for 
T = 5000 u.t. The distributions of both the MEGNO and the FT-LCN are also included. 
These figures were generated using [35j ] . 



4-1. Results at short integration times 

For the sake of making the comparison clearer, in Figs. |3] and H] we have 
changed the scale of Y through the division by the total integration time, 
T, in order to have both the MEGNO and the largest FT-LCN values of 
the same magnitude. Let us recall that a factor 2 should be added in the 
case of chaotic orbits since we are dealing with Y/T instead of Y/T (see the 
discussion at the end of Section [2]) . 

At short integration times we observe a few orbits, classified as chaotic by 
the FT-LCN, falling very close to the regular value of the MEGNO. Indeed, 
only 0.82% of the subset of chaotic orbits Oc attained MEGNO values close 
to 2.01 (~ -3.095 in Fig. El left panel) at T = 5000 u.t., while 8.82% of the 
orbits in the subset Or achieved MEGNO values within the range [2.01,10), 
indicating either their mild chaotic character or that the total integration 
time T = 5000 u.t. is not large enough for the asymptotic regular value to 
be reached. 

For the orbits in Oc, the mean of log(FT — LCN) m —1.267 and the mean 
of log(2F/T) pa —1.317, while the corresponding standard deviations are 
~ 0.617 and ~ 0.676, respectively, with a correlation coefficient of r « 0.942. 
Thus, both distributions are quite similar. 
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Figure 4: Correlations between 2Y/T and the largest FT-LCN for chaotic orbits (left 
panel) and between Y/T and the largest FT-LCN for regular orbits (right panel) for 
T = 10000 u.t. The distributions of both the MEGNO and the FT-LCN are also included. 
These figures were generated using j35|. 



For the regular sample instead, the concordance is, as expected, less for- 
tunate. The concomitant correlation coefficient r is close to 0.07. The re- 
spective mean values are FT - LCN « 0.00126 and Y/T « 0.0004, with 
rather different standard deviations, namely, 0.00016 for the distribution of 
the FT-LCNs and 0.000035 for Y/T. 

Notice must be taken of the fact that we are comparing values of the 
FT-LCNs and the MEGNO corresponding to different integration times. 
Altogether, the classification by recourse of the MEGNO provides fairly good 
results taking account that they are obtained for T = 5000 u.t., half the 
total integration time used by pj in their computation of the FT-LCN. The 
comparison of both dynamical indicators at the very same total integration 
time T is the subject of the forthcoming section. 

4-2. Results at intermediate integration times 

From Figs. H] we observe a rather good agreement between the classifica- 
tion given by the MEGNO and that due to the largest FT-LCN at T = 10000 
u.t. In fact, the chaotic component Oc appointed by the FT-LCN is re- 
attained by means of the MEGNO, i.e. all orbits in Oc have MEGNO values 
lying on the MEGNO irregularity range (above ~ —3.4 on the vertical axis 
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in the plot on the left). Since we have rescaled the MEGNO by 1/T and 2/T 
for regular and chaotic orbits respectively, those values above, but close to, 
either 0.0002 or —3.4 after 10000 u.t., suggest that the orbit could be prox- 
imate to an hyperbolic object (like unstable periodic orbits, 2D hyperbolic 
torus). On the other hand, for MEGNO values < 0.0002, the orbit may be 
close to elliptical objects (stable periodic orbits, 2D resonant elliptical torus), 
as it is shown in 

On the left panel in Fig. HJ corresponding to chaotic orbits, we distinguish 
a clearly linear correlation between log(2F/T) and log(FT — LCN) for T = 
10000 u.t. Indeed, the correlation coefficient is r « 0.95; the mean value of 
log(2F/T) w —1.29, while the concomitant mean value of log(FT — LCN) « 
-1.27. 

In the same figure on the right, we show the correlation between Y/T 
and the largest FT-LCN for the orbits in Or. In this case, the correlation 
coefficient is r w —0.026, the mean value of Y/T is rather close to 0.0002 
with a standard deviation < 10 -5 , while for the FT-LCNs the mean value is 
0.0013 with a standard deviation of order m 10~ 4 . Let us point out the sharp 
character of the distribution of the MEGNO values around the predicted one 
for regular motion, while the FT-LCNs' distribution is rather blunt, as the 
standard deviation of both distributions indicate. This fact should be deemed 
as an advantageous feature of the MEGNO over the largest FT-LCN. 

Notice must be taken that the empirical value V£ = 0.00155 adopted by 
[l[ for their classification, is greater than the mean value of the largest FT- 
LCN for regular orbits, which indeed might be a misleading factor for the 
task. 

Thus, from the regular component Or, which encompasses 1644 orbits on 
the whole, 1513 orbits have values of the Y/T in the interval [0.0001; 0.000201), 
in due accord with their stable, regular character. However, a discrepancy 
is found for a subset Or d , including the remaining 131 orbits classified as 
regular by their largest FT-LCN, whose Y/T values, however, lie within the 
range [0.000201; 0.001) revealing their possible irregular character. 

An issue to be stressed is the fact that the deemed orbits belong to differ- 
ent energy surfaces and, for each energy, a different characteristic time-scale, 
T C (E), can be defined. Moreover, in order to ensure that both the FT-LCN 
and MEGNO are well computed for a given orbit, the total integration time 
should verify T >> T C (E). Therefrom, on fixing the condition T > 10 3 T C (E) 
to obtain confident values for both indicators, we conclude that only those 
orbits with T C (E) < 10 would be properly classified for a total integration 
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Figure 5: Period of the x-axis orbit as a function of the energy adopted as T C (E) for box 
orbits, along with the approximate period of a sample of both box and tube orbits (on the 
left). Energy of the 131 orbits in Or^ (on the right). 



time of 10000 u.t. The function T C (E) for box orbits, taken as the period 
of the stable x-axis periodic orbit, is plotted in Fig. |5] on the left, which 
shows that T C (E) = 10 corresponds to an energy value E ~ —0.58. We have 
computed the approximate period of some box orbits (crosses in the figure) 
in order to test the suitability of the adopted time-scale T C (E) for boxes. 
The approximate period of some tube orbits are also included in the figure. 

The energy values corresponding to the orbits in Or^ are displayed on 
the right of Fig. El where we have labeled as group 1 those with energies 
E < —0.58, and as group 2 the ones for which E > —0.58. 

From the 131 orbits in Or^ 41 are tubes while 90 are boxes, 20 of which 
have E > —0.58, i.e. their dynamical indicators would still be in a transient 
phase. On the whole, we count 70 box orbits in Or d satisfying the condition 
T » T C {E). 

Let us remark that all orbits in Or rf have 2<F<7atT = 10000 u.t. 
so, even when they could evince some local instability, they behave as stable 
orbits from a physical point of view. Nonetheless, since our aim is to subject 
the MEGNO to a rigorous test as a dynamical indicator, we will study this 
subset of orbits in particular. 

Therefore, for the 131 orbits in Or^ we recalculate the MEGNO but for 
100000 u.t. to find that 52 orbits, having MEGNO values smaller than 3.5 at 
T = 10000, approach the regular value 2 at T = 100000, while the remaining 
79 attain greater values of the MEGNO for the larger integration time. Let 
us mention that for some orbits the MEGNO value is barely higher than 
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Figure 6: Both the largest FT-LCN provided by [l[ and the maximum LCN computed by 
our using the classical algorithm for a sample of 100 orbits in O (left). Largest FT-LCN, 
LCN, slope of the MEGNO and MEGNO for a smaller sample of 45 orbits (right). 



the adopted threshold of 2.01 (as it will be shown in section 14.4(1 . In the 
following subsection, we will identify the first group by Or^ and the second 
one by Or^. 

In order to determine the actual character of the orbits in Or^, we will 
recourse to a slight variation of the so-called Fast Lyapunov Indicator (FLI) 
(see [20]). Briefly, the FLI is defined as the supremum of the norm of the 
tangent vector 5. Thus, we will follow the evolution of (5(t)) where 

(6(t)) = - t j\(t)dt, 

with 5 = 1 1 5 1 1 and S is the solution of the variational equations, taking as 
initial value 5q = 10~ 4 . 

Since motion times larger than 10000 u.t. will be considered in the next 
section, let us first compute the LCN for a sample of 100 orbits in O for T = 



10000 u.t. using the classical algorithm of |36|. and compare the obtained 
values with the largest FT-LCN derived by pj. The result of the comparison 
is illustrated in the plot on the left of Fig. El which gives account of an 
actually quite good agreement. This encourages us to compute the LCNs 
and assimilate them with the largest FT-LCNs for T = 100000 u.t. In Fig. [6] 
on the right we also include the estimation of the largest FT-LCN derived 
from the slope of the MEGNO and the MEGNO itself for a small sample of 
orbits (45 on the whole). Let us point out that the slope of the MEGNO 
yields a better estimation of the largest FT-LCN, particularly for the regular 
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orbits, for which it lies below 10 4 , which is much smaller than both and 
Vl for T = 10000 u.t. Nonetheless, we will restrict our comparative analysis 
to the MEGNO, the LCN and the mean FLI ((£(*))), since deriving the 
expected theoretical value for the slope of the MEGNO in case of quasi- 
periodic motion is difficult. In fact, the same occurs with the mean FLI, 
for which it is not possible to determine an asymptotic value whenever the 
orbit is confined to a torus. Indeed, for such a quasi-periodic orbit, j q , the 
solution of the variational equation in M 2N can be recast as 

$ (7a(*)) ~ <5o [1 + w q (t) +t(A g + u q (t))] , 

where X q > is the linear rate of divergence around j q , and w q (t) and u q {t) 
are oscillating functions of t of bounded amplitude (in general quasi-periodic 
and with zero average), satisfying < b q < X q . The parameter X q is 

a measure of the lack of isochronicity around the orbit since it is related to 
the maximum eigenvalue of the matrix dui/dl, u and I being the frequency 
and action vectors associated to the torus, respectively (for an isochronous 
system, such as the harmonic oscillator, A = b q = for all 7). 

4-3. Results at large integration times 

In the present section we will be concerned with the temporal evolution of 
the three indicators to be compared, namely the MEGNO, the LCN and the 
mean FLI (for which no re-normalization was performed and in the case of the 
exponential growth of 8(t) the integration was stopped at (S(t)) = 10 20 ), for 
large motion times. Figs. [7] displays the typical behavior of these indicators 
for chaotic and regular orbits which are identified by their orbit number. For 
the illustration the orbits 1491, 442 and 3359 from the set Oc and 358, 2105 
and 3375 from Or — Or^ have been selected. 

The MEGNO shows a linear growth with time for the orbits in Oc, except 
for orbit 1491 for which exhibits two similar linear trends and a flat behavior 
between ~ 20000 and ~ 60000 u.t., suggesting that during this time interval 
the orbit may be close to some elliptic structure. On the other hand, for the 
orbits in Or — Or^ the MEGNO asymptotically approaches the predicted 
value, 2, both at 10000 u.t and 100000 u.t. 

In regards to the LCN, a similar behavior is observed for the orbit 1491, 
while for the regular orbits converges to the theoretical expected value V* = 
InT/T ~ 1.2 x 10~ 4 . 

Finally, for those orbits in Oc the mean FLI displays an almost expo- 
nential dependence with time (in fact it goes as e at /t), while it attains much 
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Figure 7: Typical temporal evolution of the MEGNO, the LCN and the mean FLI for 
chaotic and regular orbits (e.g. in Oc and Or — Or^, respectively). For the plots in the 
last panel, the runs were stopped at (S(t)} = 10 20 in order to avoid overflow. 



lower values for orbits in Or — Or^ for which it depends with time in a linear 
fashion (note the logarithmic scale in the vertical axis), 

Thus, these figures provide information about the expected behavior of 
the three indicators in the cases of both regular and chaotic motion, which 
will be of use to determine the character of those orbits in Or^ and Or^. 
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Let us recall that we will restrict our study to those orbits satisfying the 
condition E < —0.58 and for which a good estimation of their period is at 
hand. 

The correlation between the values of MEGNO and the LCN at T = 
100000 u.t. for orbits in Or^ and Or^ is presented in Fig. [HJ For the 
orbits in Or^, the mean value of log(2Y/T) ~ —3.28 while the mean of 
log(LCN) ~ —3.13, and the standard deviations are 0.72 and 0.56 respec- 
tively, the correlation coefficient being close to 0.98. Therefore, not only the 
correlation between both indicators is quite good, but the two first moments 
of their concomitant distributions are rather similar as well. Let us mention 
that, though we are computing the MEGNO for very large times, most of the 
orbits in Or^ attain values in the range —4.5 < log(2F/T) < —2.5, reveal- 
ing that these orbits are mild chaotic, for which the mean Lyapunov time is 
Tiyap ~ 1500 u.t., their characteristic period being smaller than 10 u.t. 



t-2.5 



-3 -2.5 -2 -1.5 

log[LCN] 



1 .98e-05 
t 1.976-05 



i Orbit 21 22 ■ 

0.00012 0.00014 0.00016 0.00018 
LCN 



Figure 8: Correlations between 2Y/T and the LCN for orbits in Or^ (left panel) and 
between Y/T and the LCN for orbits in Or s d (right panel) at T = 100000 u.t. The orbits 
marked with a full dot will be studied separately, following the temporal evolution of the 
three dynamical indicators. 



Meanwhile, for orbits in Or^, the mean value of Y/T m 2 x 10 -5 with 
a standard deviation of 1.3 x 10~ 7 , while the mean LCN is close to 0.00014 
with a standard deviation of about 1.7 x 10~ 5 . The corresponding correlation 
coefficient is 0.1. Again we point out the sharp distribution of Y/T around 
the expected theoretical value. The standard deviation of both distributions 
differ in two orders of magnitude (notice should be taken of the different 
scales onto the vertical and horizontal axis in the right plot of Fig. |8]). At T = 
100000 u.t. we attain values of Y/T that provides a fairly good estimation 
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of the true LCN, a = 0, namely, of order 10 -5 , while the LCN computed by 
recourse to the standard algorithm barely abuts 10 -4 . 

The explanation is clear; the factor InT w 11 in V* is the responsible for 
this slower convergence of <j\ to a = as T — > oo, and Y/T tends to o faster 
than In T/T. Indeed, for stable motion and integration times of order of, or 
larger than T « 22000 u.t., \Y/T\ < 10~ 4 , while the LCN computed using 
the standard algorithm yields |o~i| < 10~ 3 . In fact, 

Y/T 2 
o"i InT 

Finally, the full dots in Fig. [8] correspond to five orbits selected as sam- 
ples of Or^ and Or^, for which the study of the temporal evolution of the 
three indicators for T = 100000 u.t. will serve to determine their dynamical 
behavior. This issue will be undertaken in the forthcoming subsection. 

4-4- On the Ord orbits 

Let us be concerned with the detailed study of some sample orbits of 
Or^. First we will aim our attention at orbits 2881, 537 and 3297 belonging 
to Or^, and 2122 G Or^. The study of orbit 228 G Or^ will be addressed 
separately. 

Fig. [9] displays the temporal evolution of the three indicators for both 
integration times, namely, T = 10000 and T = 100000 u.t., corresponding 
to the selected sample orbits. It can clearly be observed that for T = 10000 
u.t. almost all orbits exhibit a stable behavior. Yet, both the MEGNO and 
the mean FLI evince an incipient increase for orbits 537 and 3297, which is 
missed by the LCN. As time increases, both orbits clearly separate from the 
rest, all the three indicators giving account of this fact. 

Meanwhile, orbits 2122 G Or^ and 2881 G Or^ seem to evolve in a similar 
fashion. Though, the final MEGNO value for orbit 2122 is sharply 2, while 
for 2881 is slightly above the regular value, which might indicate a rather 
mild unstable character of this orbit (see discussion below). 

In Fig [10] we present the time evolution of both the MEGNO and the 
LCN, on the left and right panel respectively, for orbit 228 G Or^. It is 
interesting to note the particular behavior of the indicators for this orbit. 
From the plot on the left, the trajectory looks like a stable quasi-periodic 
orbit up to t < 2000 u.t., then the MEGNO grows linearly for a rather short 
time interval to reach a nearly constant value, around 2.6, and attains the 
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Figure 9: Time evolution of the MEGNO, the LCN and the mean FLI for the sample 
orbits in Or^ and Or^ for T = 10000 and T = 100000 u.t. In the bottom panels we 
separately plot the evolution of the mean FLI for different times. 




Figure 10: Time evolution of the MEGNO and the LCN for orbit 228 e Or^ for T = 10000 
and T = 100000 u.t. 



value 3 at T = 10000 u.t. Note that the LCN at this time is very close to 
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the theoretical expected one, around 0.001. When the integration time is 
increased the MEGNO grows up to higher values but not in a linear fashion, 
while the LCN seems to decrease, though it approaches a larger value than 
the one corresponding to regular motion at T = 100000 u.t. The peculiar 
behavior of the indicators for this orbit encourages a more detailed study of 
its neigbourhood in phase space in order to grasp its actual dynamical nature. 
This chore will be performed by analysing its immediate neihgbourhood in 
phase space, which might provide us with valuable dynamical information. 
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Figure 11: MEGNO and the LCN for T = 10000 and T = 100000 u.t. for 125 orbits in 
a domain of size 10~ 7 around orbit 228. The last value corresponds to the mean value of 
the indicators, also represented by the horizontal line. 



Thus, in Fig. [TJJwe present the MEGNO and LCN values at T = 10000 
and T = 100000 u.t. for a set of 125 orbits taken at random in a neighborhood 
of size 10~ 7 centered at orbit 228, whose concomitant values are depicted by 
full dots in each plot. Note that at T = 10000 the MEGNO for this orbit is 
very close to 2, while the mean value of the indicator for this set of orbits 
is about 11. On increasing the integration time, it becomes quite clear that 
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orbit 228 is in fact chaotic. A similar behavior is observed in regards to the 
LCN values. Altogether, the figure suggests that this orbit might lie in a 
complex dynamical region of phase space. 

In order to confirm this conjecture, in Fig.[12]we show a MEGNO contour 
plot in the momenta space for the exact energy value of orbit 228, constructed 
by taking as initial conditions the position of the orbit in configuration space 
and (p x ,Pz) varying over a grid of about 10 6 points. The MEGNO values 
correspond to T = 1000 u.t., which turns out to be a proper final time of 
integration since T « 10 3 T C (E), and the characteristic time scale is about 1 
for an energy E ~ —3 (as follows from Fig. [5]). 

The dark regions corresponds to strong chaotic motion for which Y > 20, 
while the white ones with Y < 2.01 reveal stable motion. Light gray zones 
refer to slight unstable (or even regular) motion, 2.01 < Y < 3 and the dark 
gray regions represent mild chaotic motion 3 < Y < 20. This plot reveals the 
complex resonance structure of phase space when projected onto the plane 
ipxiPz) a t this energy level, and we can clearly see that orbit 228 lies inside 
a resonance crossing. This should explain the pathological behavior of this 
orbit and its surroundings. 

Just to end this section, let us investigate the surroundings of orbit 2881 
which belongs to Or^. For that sake, let us consider 125 orbits selected at 
random within a neibourhood of size 10~ 7 and compute both their MEGNO 
and LCN at T = 10000 and T = 100000 u.t. The results are displayed in 
Fig. O to show that the MEGNO values lie in the range 1.95 < Y < 3 for 
T = 10000, and on considering larger motion times the MEGNO interval gets 
even narrower, e.g. (1.98,2.27) for T = 100000 u.t. Therefore, though the 
MEGNO for orbit 2881 is slightly higher than the threshold 2.01, namely, 
2.014 at T = 100000, the orbit should be considered at all means stable. 
An analogous result provides the LCN, that has been computed for the very 
same orbits, including orbit 2881 for which, at T = 10000 u.t., attains a value 
less than the empirical critical value V£ adopted by jlj, but higher than V*. 

Finally, it would be interesting to consider the MEGNO values of all 
orbits in Or^ for T = 100000 u.t., which are displayed in Fig. HU It can 
there be noticed that at least 14 orbits in Or^ should actually be included 
in the set Or^ (e.g. 2881), since their MEGNO values are rather too close 
to the regular value 2. In fact, it might be inaccurate for them to consider 
the factor 2 in 2Y/T, necessary in the case of chaotic orbits, on looking for 
correlations with the LCN, since these orbits do not increase linearly with 
time. 
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" 0.1 0.2 0.3 0.4 0.5 0.6 0.7 O.B 

Figure 12: MEGNO contour plot for a grid of 10 6 initial conditions in momenta space 
(p x ,Pz) for E = —3.041052762 and initial position of orbit 228. The total motion time 
is T — 1000 « 10 3 T C (E). Black (Y > 20) corresponds to strong chaotic zones, white 
indicates regular regions Y < 2.01, light gray to slight unstable or even some regular 
orbits (2.01 < Y < 3) and dark gray indicates mild chaotic domains (3 < Y < 20). The 
arrow indicates the location of orbit 228. 



On transfering these 14 orbits from Or^ to Or^, the recomputation of 
the concomitant resulting distributions deliver for Or^, the mean value of 
log(2y/T) « —3.04 with a standard deviation of about 0.55 and the mean 
of log(LCN) PS —2.98 with a standard deviation close to 0.51, the correlation 
coefficient being r ps 0.99. For Or^ there results a mean value of Y/T ~ 
2 x 10~ 5 with a standard deviation of ps 2.3 x 10~ 7 , while the mean LCN 
is ~ 1.4 x 10~ 4 with a standard deviation of ps 1.8 x 10~ 5 . The correlation 
coefficient in this case barely amounts 0.22, indicating again no correlation 
between both indicators for regular orbits. Let us notice that no significant 
changes arise as a consequence of the transposition performed. 

5. Discussion 

We have shown that the MEGNO is a suitable fast indicator to separate 
regular from chaotic motion. Further, it is particularly useful to investigate 
the nature of orbits that have a small but positive Lyapunov number. 

Besides we have shown a rather good correlation between the MEGNO 
and the FT-LCN values for short, moderate and large integration times for 
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Figure 13: MEGNO and the LCN for T = 10000 and T = 100000 u.t. for 125 orbits in a 
domain of size 10~ 7 around orbit 2881. The last value corresponds to the mean value of 
the indicators, also represented by the horizontal line. 
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Figure 14: MEGNO values for the 79 orbits in Or^ for T = 100000 (left). Zoom around 
Y = 2 (right). 



chaotic orbits, while the MEGNO provides better results for regular motion. 
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In fact, it has the advantage that Y/T converges to the null value of a 
faster than the classical algorithm to compute the LCN. Another recourse to 
derive low values for the LCN in the case of quasi-periodic motion consists 
in computing the slope of the MEGNO. 

The FLI looks also as a reliable fast indicator, but it does not provide 
any reference value for regular motion, so it may be useful to explore phase 
space rather than to investigate the nature of a given orbit, unless we follow 
the time evolution of this indicator. 

Finally, though Y — > 2 when t — > oo and for quasi-periodic motion its 
convergence would be very fast, one should keep in mind that a single orbit 
is in general dominated by the dynamics of its surroundings. Therefore if 
it is necessary to determine the strict character of that orbit, the threshold 
value of a regular orbit would be Y ~ d where d could be taken as the mean 
value of Y over a rather small domain around the orbit. 
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